#ifndef UHM_MAT_UHM_HXX_
#define UHM_MAT_UHM_HXX_

#define UHM_MAT_UHM_UNROLL 4

namespace uhm {

  class UHM_Base_ {
  public:
    enum { KIND_A_BEGIN = 1,
           KIND_ATL = 1,
           KIND_ATR, 
           KIND_ABL,
           KIND_ABR,
           KIND_A_END,
           KIND_B_BEGIN = 11,
           KIND_XT = 11,
           KIND_XB,
           KIND_BT,
           KIND_BB,
           KIND_RT,
           KIND_RB,
           KIND_B_END,
           KIND_P = 101,
           KIND_QR = 201} Mat_Kind_;
  private:
  protected:
    Int_ _rhs, _blk;

  public:
    UHM_Base_() : _rhs(0) { }
    UHM_Base_(Int_ rhs, Int_ blk) : _rhs(rhs), _blk(blk) { }
    UHM_Base_(const UHM_Base_ &b) : _rhs(b._rhs), _blk(b._blk) { }

    virtual ~UHM_Base_() { }
    
    UHM_Base_& operator=(const UHM_Base_ &b) {
      this->_rhs = b._rhs;
      this->_blk = b._blk;

      return (*this);
    }

    void set(Int_ rhs, Int_ blk) {
      _rhs = rhs; _blk = blk;
    }

    Int_ get_rhs() const { return _rhs; }
    Int_ get_blk() const { return _blk; }

    virtual Double_ get_leading_flop() { return 0.0; }
  };
  
  template<class A_, class B_>
  class UHM_ : public Mat_Base_,
               public UHM_Base_, 
               public Disp_ {
  private:
  protected:
    A_ _ATL, _ATR, _ABL, _ABR, _P;
    B_ _XT, _XB, _BT, _BB, _RT, _RB;

    //Bool_ _pivotable;

  public:
    typedef A_ _A_;
    typedef B_ _B_;

    UHM_() : 
      Mat_Base_(),
      UHM_Base_() { }

    UHM_(Int_ type, Int_ uplo, Int_ sym, 
         Int_ mt, Int_ mb, 
         Int_ rhs, Int_ blk) :
      Mat_Base_(type, uplo, sym, mt, mb), 
      UHM_Base_(rhs, blk),
      _ATL(type, uplo,     sym,       mt, mt),
      _ATR(type, FULL, UNSYM, mt, mb),
      _ABL(type, FULL, UNSYM, mb, mt),
      _ABR(type, uplo,     sym,       mb, mb),
      _XT (type, FULL, UNSYM, mt, rhs), 
      _XB (type, FULL, UNSYM, mb, rhs), 
      _BT (type, FULL, UNSYM, mt, rhs), 
      _BB (type, FULL, UNSYM, mb, rhs), 
      _RT (type, FULL, UNSYM, mt, rhs), 
      _RB (type, FULL, UNSYM, mb, rhs),
      _P  (INT,  FULL, UNSYM, mt, 1) { }
    
    UHM_(Int_ type, Int_ mt, Int_ mb, Int_ rhs, Int_ blk) :
      Mat_Base_(type, mt, mb), 
      UHM_Base_(rhs, blk),
      _ATL(type, FULL, UNSYM, mt, mt),
      _ATR(type, FULL, UNSYM, mt, mb),
      _ABL(type, FULL, UNSYM, mb, mt),
      _ABR(type, FULL, UNSYM, mb, mb),
      _XT (type, FULL, UNSYM, mt, rhs),
      _XB (type, FULL, UNSYM, mb, rhs),
      _BT (type, FULL, UNSYM, mt, rhs),
      _BB (type, FULL, UNSYM, mb, rhs),
      _RT (type, FULL, UNSYM, mt, rhs),
      _RB (type, FULL, UNSYM, mb, rhs),
      _P  (INT,  FULL, UNSYM, mt, 1) { }

    UHM_(const UHM_ &b) : 
      Mat_Base_(b), 
      UHM_Base_(b),
      _ATL(b._ATL),     _ATR(b._ATR),
      _ABL(b._ABL),     _ABR(b._ABR),
      _XT (b._XT),      _XB (b._XB),
      _BT (b._BT),      _BB (b._BB),
      _RT (b._RT),      _RB (b._RB),
      _P  (b._P) { }

    virtual ~UHM_() { free_all(); }

    UHM_& operator=(const UHM_ &b) {

      Mat_Base_::operator=(b);
      UHM_Base_::operator=(b);

      this->_ATL = b._ATL;      this->_ATR = b._ATR;
      this->_ABL = b._ABL;      this->_ABR = b._ABR;
      
      this->_XT = b._XT; this->_BT = b._BT; this->_RT = b._RT;
      this->_XB = b._XB; this->_BB = b._BB; this->_RB = b._RB;

      this->_P  = b._P;

      return (*this);
    }


    virtual void set(Int_ type, Int_ uplo, Int_ sym,
                     Int_ mt, Int_ mb,
                     Int_ rhs, Int_ blk) {
      this->Mat_Base_::set(type, uplo, sym, mt, mb);
      this->UHM_Base_::set(rhs, blk);
      _ATL.set(type, uplo, sym,   mt, mt);
      _ATR.set(type, FULL, UNSYM, mt, mb);
      _ABL.set(type, FULL, UNSYM, mb, mt);
      _ABR.set(type, uplo, sym,   mb, mb);
      _XT.set (type, FULL, UNSYM, mt, rhs);
      _XB.set (type, FULL, UNSYM, mb, rhs);
      _BT.set (type, FULL, UNSYM, mt, rhs);
      _BB.set (type, FULL, UNSYM, mb, rhs);
      _RT.set (type, FULL, UNSYM, mt, rhs);
      _RB.set (type, FULL, UNSYM, mb, rhs); 

      _P.set (INT, FULL, UNSYM, mt, 1); 
    }

    virtual void set(Int_ type, Int_ mt, Int_ mb, Int_ rhs, Int_ blk) {
      this->set(type, FULL, UNSYM, mt, mb, rhs, blk);
    }


    Int_ get_mt() const { return this->_m; }
    Int_ get_mb() const { return this->_n; }

    A_& ATL() { return _ATL; }
    A_& ATR() { return _ATR; }
    A_& ABL() { return _ABL; }
    A_& ABR() { return _ABR; }

    B_& XT() { return _XT; }
    B_& XB() { return _XB; }
    B_& BT() { return _BT; }
    B_& BB() { return _BB; }
    B_& RT() { return _RT; }
    B_& RB() { return _RB; }

    A_& P() { return _P; }

    virtual A_& A(Int_ kind) {
      switch (kind) {
      case KIND_ATL: return _ATL;
      case KIND_ATR: return _ATR;
      case KIND_ABL: return _ABL;
      case KIND_ABR: return _ABR;
      }

      std::cout << "Given Kind A = " << kind << "\n";
      UHM_ERROR(false, UHM_ERROR_MAT_KIND);

      // dummy return...
      return _ATL;
    }

    virtual B_& B(Int_ kind) {
      switch (kind) {
      case KIND_XT: return _XT;
      case KIND_XB: return _XB;
      case KIND_BT: return _BT;
      case KIND_BB: return _BB;
      case KIND_RT: return _RT;
      case KIND_RB: return _RB;
      }

      std::cout << "Given Kind B = " << kind << "\n";
      UHM_ERROR(false, UHM_ERROR_MAT_KIND);

      // dummy return...
      return _XT;
    }

    virtual Double_ get_storage()  {
      Double_ r_val=0.0;
#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_B_BEGIN;i<KIND_B_END;++i)
        r_val += B(i).get_storage();

#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_A_BEGIN;i<KIND_A_END;++i)
        r_val += A(i).get_storage();

      return r_val;
    }

    virtual Double_ get_allocated_storage()  {
      Double_ r_val=0.0;
#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_B_BEGIN;i<KIND_B_END;++i)
        r_val += (B(i).get_storage() * !B(i).is_flat_buffer_null());

#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_A_BEGIN;i<KIND_A_END;++i)
        r_val += (A(i).get_storage() * !A(i).is_flat_buffer_null());

      return r_val;
    }

    virtual void create_flat_without_buffer() { 
#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_B_BEGIN;i<KIND_B_END;++i) 
        B(i).create_flat_without_buffer();

#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_A_BEGIN;i<KIND_A_END;++i) 
        A(i).create_flat_without_buffer();
      
      _P.create_flat_without_buffer();
    }

    virtual void create_flat_buffer() {
      create_flat_buffer(0);
    }
    virtual void create_flat_buffer(Int_ abr) {

#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_B_BEGIN;i<KIND_B_END;++i) 
        B(i).create_flat_buffer();

      _ATL.create_flat_buffer();
      if (abr)
        _ABR.create_flat_buffer();

      switch(this->_uplo) {
      case FULL:
        _ATR.create_flat_buffer();
        _ABL.create_flat_buffer();
        break;
      case LOWER_TRIANGULAR:
        _ABL.create_flat_buffer();
        break;
      case UPPER_TRIANGULAR:
        _ATR.create_flat_buffer();
        break;
      }

      _P.create_flat_buffer();
    }

    virtual void create_flat_ABR_buffer() {
      _ABR.create_flat_buffer();
    }

    virtual void free_flat_ABR_buffer() {
      _ABR.free_flat_buffer();
    }

    virtual void free_flat_buffer() {
#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_B_BEGIN;i<KIND_B_END;++i) 
        B(i).free_flat_buffer();

#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_A_BEGIN;i<KIND_A_END;++i)
        if (!A(i).is_flat_buffer_null())
          A(i).free_flat_buffer();

      _P.free_flat_buffer();
    }

    virtual void create_hier() {
#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_B_BEGIN;i<KIND_B_END;++i) 
        B(i).create_hier(_blk, _blk);

#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_A_BEGIN;i<KIND_A_END;++i) 
        A(i).create_hier(_blk, _blk);

      _P.create_hier(_blk, _blk);
    }

    virtual void free_hier() {

#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_B_BEGIN;i<KIND_B_END;++i) 
        B(i).free_hier();

#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_A_BEGIN;i<KIND_A_END;++i)
        A(i).free_hier();

      _P.free_hier();
    }

    virtual void free_all() {
#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_B_BEGIN;i<KIND_B_END;++i) 
        B(i).free_all();

#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_A_BEGIN;i<KIND_A_END;++i)
        A(i).free_all();

      _P.free_all();
    }

    virtual void set_zero_A() {
#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_A_BEGIN;i<KIND_A_END;++i)
        A(i).set_zero();

      _P.set_zero();
    }

    virtual void set_zero_B() {
#pragma unroll(UHM_MAT_UHM_UNROLL)
      for (Int_ i=KIND_B_BEGIN;i<KIND_B_END;++i)
        B(i).set_zero();
    }

    virtual void set_zero() {
      this->set_zero_A(); this->set_zero_B();
    }
    virtual void set_identity() {
      _ATL.set_identity();
      _ATR.set_zero();
      _ABL.set_zero();
      _ABR.set_identity();
    }
    virtual std::ostream& content_general(std::ostream& os) const {
      os << "- UHM General Info, " 
         << " Type " << this->_type  
         << " Uplo " << this->_uplo
         << " Sym " << this->_sym
         << std::endl;
      return os;
    }

    virtual std::ostream& content_A(std::ostream& os) const  {
      os << "- ATL - " << std::endl;
      os << _ATL.get_flat() << std::endl;

      os << "- ATR - " << std::endl;
      os << _ATR.get_flat() << std::endl;

      os << "- ABL - " << std::endl;
      os << _ABL.get_flat() << std::endl;

      os << "- ABR - " << std::endl;
      os << _ABR.get_flat() << std::endl;

      os << "- P - " << std::endl;
      os << _P.get_flat() << std::endl;

      return os;
    }

    virtual std::ostream& content_B(std::ostream& os) const  {
      os << "- XT -" << std::endl;
      os << _XT.get_flat() << std::endl;

      os << "- XB -" << std::endl;
      os << _XB.get_flat() << std::endl;

      os << "- BT -" << std::endl;
      os << _BT.get_flat() << std::endl;

      os << "- BB -" << std::endl;
      os << _BB.get_flat() << std::endl;

      os << "- RT - " << std::endl;
      os << _RT.get_flat() << std::endl;

      os << "- RB - " << std::endl;
      os << _RB.get_flat() << std::endl;

      return os;
    }

    virtual std::ostream& content(std::ostream& os) const  {
      this->Disp_::content(os);
      content_general(os);
      content_A(os);
      content_B(os);

      return (os << std::endl);
    }
  };
}

#endif
